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ABSTRACT 



The longitudinal stability of an aircraft at or near 
stall was examined using the digital computer as an experi- 
mental tool to solve the longitudinal equations of motion. 

A linear analysis determined the effect of lift curve slope 
variation. An investigation was made to identify the non- 
linear lift curve variations needed to create the often 
observed "rocking-chair" or "porpoising" stall trait. The 
characteristics of this limit-cycle oscillation were examined. 
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I. 



INTRODUCTION 



Federal and military aviation requirements stipulate that 
all aircraft possess stalling characteristics which comply 
with established criteria. The Federal Aviation Regulations, 
Part 23 of "Airworthiness Standards" state, in part, that 
"acceptable stalling characteristics be demonstrated. .. in 
straight flight with wings level .. .where the primary control 
manipulation is a steady progressive upward movement of the 
elevator until the aircraft is stalled." A desirable stall 
characteristic in a longitudinal sense would include a small 
but distinguishable drop in the nose of the aircraft, which 
can not be overridden by the pilot. Also, if the elevator 
control were eased forward, the aircraft should promptly 
return to unstalled flight. 

The exact definition of aircraft stall will differ among 
the pilot, wing design specialist, stability and control ex- 
pert, and all others who will have their own interpretation. 
Highly sweptwinged, supersonic aircraft have discredited the 
idea of flow separation, Ct i and a stable nose-down 

pitching moment all occur ing simultaneously. The Federal 

\ 

Aviation Regulations, Part 25, defines the calibrated stalling 
speed for a civil air transport as "...the minimum steady speed, 
in knots, at which the aircraft is« controllable ... " For the 
purpose of this treatise, stall may be described as the con- 
dition when flow breakdown, primarily over the main wing, 
causes significant nonlinear effects in the moment and lift 
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characteristics. The effect of an asymmetric flow condition 
leading to wing drop will not be considered, since only the 
longitudinal stability traits are being analyzed. 

A considerable amount of literature is available on the 
subject of aircraft static stability at stall, and it de- 
scribes the aircraft's initial tendencies, or stall charac- 
teristics. Some literature is available concerning dynamic 
stability at deep stall, or large angles of attack. There 
is, however, very little treatment of airplane dynamics in 
the area of , of post-stall time histories in that 

narrow region of the stall break, or of the oscillations 
that might occur. 

Most pilots with any degree of proficiency have ex- 
perienced a "rocking-chair" or "porpoising" type of stall 
where, with constant elevator control, the aircraft's nose 
oscillates vertically about some point as the wing stalls, 
recovers, rotates back up into the stall, and so on. It 
was the purpose of this study, therefore, to numerically 
investigate the mechanism required to create this type of 
limit-cycle, or possibly divergent, rocking-chair motion, 
using the complete set of aircraft equations of motion, and 
to attempt to describe the major factors influencing this 
phenomenon . 

Because of time considerations and the availability of 
data, this study was limited to only the longitudinal sta- 
bility, with no roll, yaw, or coupling considered, of a 
small, straight-winged jet aircraft with zero thrust. All 
computations were accomplished at the W.R. Church Computer 
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Center, Naval Postgraduate School, Monterey, California, 
using the IBM 360/67 digital system and standard FORTRAN 
language. The study was accomplished during the period 
September 1971 through June 1972. 

II. AIRCRAFT MODEL AND ASSUMPTIONS 

The input data used in this study were those for an F-94 
Aircraft obtained from Ref. 1. The aircraft is a straight- 
winged, tandem-seated, single-engine jet, with the center of 
gravity positioned well within the fore and aft limits. A 
zero thrust landing configuration at sea level on a standard 
day exists throughout this analysis. The aircraft has 
basically linear aerodynamic characteristics except for lift, 
drag and pitching moment coefficients, which are all functions 
of angle of attack. Terms relating to such variables as 
elevator control effectiveness and pitch and angle of attack 
damping are also assumed to be linear. The aircraft was 
considered stick-fixed, with the elevator being the only means 
of longitudinal control and trim. The case of pitch-up at 
stall was neglected. Figure 1 shows the aircraft axis system 
and the related variables. 

III. DISCUSSION AND RESULTS 

Static longitudinal stability may be considered as the 
tendency of the aircraft to return to static equilibrium, 
while dynamic longitudinal stability considers the resulting 
aircraft motion as a function of time. The existence of 
static stability does not necessarily intimate the existence 
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of dynamic stability. However, the existence of asymptotic 
dynamic stability in the sense of Lyapunov (Ref. 4) implies 
static stability. 

The factors influencing aircraft stability at a stalled 
condition must be known in order to predict or create an 
instability. The airplane has six degrees of freedom: ro- 
tation in roll, pitch and yaw, and translation in the 
horizontal, vertical and lateral directions. The principal 
variables in the longitudinal motion are: 

(1) the pitch attitude of the airplane, 

(2) the angle of attack, 

(3) the flight velocity. 

The equations of motion describing the longitudinal motion of 

the aircraft are found in Ref. 2. These equations are or- 
dinary nonlinear differential equations in five variables, 
with time as the independent variable. A complete derivation 
of the equations and their application to the problem of 
aircraft stall may be found in Ref. 3. 

To perform a sensitivity analysis on the stability de- 
rivatives near stall, the equations were simplified by as- 
suming a constant elevator position and a horizontal initial 
flight path. Reference 2, Chapter 6, reduces the equations 
to three homogeneous algebraic equations in the unknowns 
u, o(, and B. By setting the determinant of the coeffi- 
cients equal to zero, the roots, A r of the characteristic 
equation may be found. Expansion of this "stability deter- 
minant" results in a quartic equation for A of the form: 
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Ax'" +Bx^ *C?^ *D\ +E=0 



The solution of this equation results in two pairs of complex 
conjugate roots. These roots may be expressed either in non- 
dimensional time or real time, with our interest focusing 
on the damping ratio, and the natural frequency. 

These numbers are easily calculated along with the resulting 
period of oscillation. For a dynamically stable aircraft, 
the real part of all of the roots must be negative; that is, 
the damping ratio is positive. An instability results if 
any of the daimping terms are negative. 

Upon analysis of the equations of motion, the only sta- 
bility derivatives making a significant contribution to 
stability were: , and . By varying these 

parameters and observing the damping ratios of the resulting 
roots of the stability determinant, a feel for the important 
factors in aircraft stability may be developed. 

Only two stability derivatives are of sufficient magnitude 
to contribute to dynamic stability. The first of these is 
C , the static stability derivative. As C becomes 
positive, the aircraft becomes statically unstable. This 
condition will not be considered since most aircraft demon- 
strate positive static stability through the stall, and dy- 
namic stability is of concern here. The second stability 
derivative of concern is , the lift curve slope. The 

results of varying lift curve slope are shown in Table 1. 
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The aircraft demonstrates dynamic stability when is 

positive, zero or slightly negative. When the lift curve 
slope reaches approximately -5.0, one pair of the roots of 
the characteristic equation has a positive real part which 
indicates aircraft instability. The instability remains as 
the lift curve slope is increased negatively (c.f., Fig. 2), 

A negative slope of the lift curve could occur after the 
maximum lift is reached in the area of stall, or in the case 
of a sharp stall break occurring locally. 

It must be noted at this point that only the short period 
motion was affected by this analysis. One would expect to 
find, therefore, that the aircraft has a stable long period 
or phugoid mode, but is unstable in angle of attack, pitch 
angle and normal acceleration (g's), with the velocity re- 
maining relatively constant as is characteristic of the short 
period mode. 

The computer program used to mathematically "fly" the 
aircraft is discussed at length in Ref. 3. Basically, there 
are four variables affecting the longitudinal motion: the 
resultant aircraft velocity, V, the pitch angle, 0, the rate 
of change of pitch angle, 6 , and the angle of attack, 0< . 

All are functions of time, t. The rate of change of angle 
of attack is considered small and is neglected. Five dif- 
ferential equations were used to describe the longitudinal 
motion, the fifth occurring by resolving the resultant air- 
craft velocity, V, into its components, u and w, in the x and 
z directions, respectively. The equations were solved by a 
fourth-order Runge-Kutta integration scheme. Input data' 
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consist of the aircraft parameters such as wing area, moments 
of inertia, thrust coefficient, air density, etc., and a 
table of data with aircraft lift, drag and pitching moment 
coefficients as functions of angle of attack. These input 
values remain constant throughout any maneuvering the air- 
craft is required to perform. The only initial conditions 
read into the program are the aircraft gross weight, desired 
angle of attack and angle of pitch at trim. 

Subroutine. SPLINl calculates the value of Ct , C and C 

L D m 

for the angle of attack specified from the data table by 
interpolation using a cubic curve-fitting scheme. Subroutine 
TRIM uses these data for a calculation of the corresponding 
velocities, elevator angle and thrust using the force and 
moment equations at static equilibrium. These calculated 
data are subsequently used as the initial conditions to start 
the Runge-Kutta integration scheme. Aircraft maneuvering 
and trim is performed by logic statements controlling elevator 
angle. Appendix C presents the entire computer program. 

The effect of a large negative value of lift curve slope 
(C =-10.0 per radian) beyond the stall was numerically 
introduced in the table look-up data for angles of attack 
greater than 0.4084 radians. Angles of attack less than this 
corresponded to the basic airplane situation as shown in 
Figure 3a. ' For angles of attack greater than this value, the 
lift curve slope was linear with a negative slope of -10.0 
per radian. This type of drop-off in lift coefficient is not 
unlike that encountered in wing leading edge flow separation 
or stall. 
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The aircraft was flown (numerically on the computer) at 
an initial trim condition of cX. =0.4083 radians (very near 
C ), 0=0.0989 radians, and a gross weight of 12,359 lbs. 

for three seconds. Then a trailing-edge-up elevator impulse 
(back stick) of 0.50 degrees (0.008725 radians) was inserted 
for 0,05 seconds. The resulting time histories show the 
aircraft to exhibit divergent oscillations, as predicted 
in Fig 2, in angle of attack, pitch angle, and lift coef- 
ficient. The aircraft velocity also shows divergence, al- 
though the amplitude was very small compared with the other 
parameters, thus confirming the short period tendencies. 

The case of stall with a smooth lift curve as defined 
by Fig. 3a was investigated. Aircraft time histories with 
the same initial conditions are shown in Figures 3b through 
3e. This type of lift curve is typical of a wing trailing 
edge flow separation or stall. The maximum negative slope 
reached was =-2.0. The damping ratio for this case 

was positive, with negative real parts for the roots of the 
characteristic equation, indicating an asymptotic type of 
dynamic stability in the sense of Lyapunov, with a dying out 
of the oscillatory motion to a static equilibrium condition. 
Mathematically, the aircraft model and the program used to 
manipulate it have demonstrated both dynamic stability and 
instability in the region of stall based on the slope of the 
lift curve at angles of attack beyond that corresponding to 
maximum lift. 

The problem of "rocking-chair" stall would indicate a 
limit-cycle type of oscillation. Reference 4 describes a 
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limit cycle as a nonlinear oscillation ... "due to the presence 
of nonlinear terms in the differential equation." A stable 
limit cycle is characterized as being a unique closed curve 
in the phase plane, to which all other nonclosed trajectories 
approach in the form of spirals winding onto the limit cycle. 

It is a self-sustained oscillation, independent of the initial 
conditions. The differential equations governing the longi- 
tudinal motion are nonlinear as well as the lift, drag and 
moment coefficients being functions of angle of attack. With 
an instability in the system at static equilibrium due to the 
negative lift curve slope, it is possible for a limit cycle 
to exist. 

To prevent an oscillatory divergence as shown previously, 
the instability was bounded by a stable lift curve slope. It 
may be likened to a reverse relay in a nonlinear control 
system. For simplicity in programming, a vertical (or infinite) 
lift curve slope was used locally. This allows the original 
lift curve table look-up data to be used and a simple logic 
statement decreases the interpolated value of lift coefficient 
by an amount ACj^ if the angle of attack were greater than a 
specified amount. Appendix D shows the function subprogram 
FI and the subroutine TRIM modified from Appendix C, each by 
only one statement. Figure 4 shows the resulting lift curve 
generated by this modification. 

Using the same initial conditions from the previous two 
runs, the same 1/2-degree elevator impulse at three seconds, 
and changing only the lift curve to that in Fig. 4, a limit 
cycle oscillation, or "rocking-chair" stall was generated-. 



17 



The AC used was 0.050. Figures 5b through 5e show lift 

Xj 

coefficient, angle of attack, pitch angle and velocity as 
functions of time. The small elevator impulse begins an 
oscillation which grows to a maximum amplitude for the dura- 
tion of the flight time. Figure 5a shows the phase plane 
plot of pitch angle ( © ) versus pitch rate ( © ) . The curve 
winds onto a closed path trajectory, indicative of a limit 
cycle . 

Amplitude of the limit cycle oxcillation could be control- 
led by the magnitude of the AC on the lift curve. If 

Xi 

^C were reduced by one-half, from 0.050 to 0.025, the 
peak amplitude of the oscillation would be halved. Conversely, 
if were doubled to 0.10, the amplitude of the limit 

cycle would be doubled. The period was not noticeably af- 
fected. Figure 6 shows the time history of the pitch angle 
for varied A Cj^ values. The initial conditions and elevator 
impulse remains the same as for previous runs. 

The magnitude of the elevator impulse had little effect 
on the limit cycle oscillation. Figure 7 shows the time 
histories of pitch angle for impulses of 0.05 degrees and 2.0 
degrees. The amplitude and period of the limit cycle remain 
the same.-' The only change was within the first fifteen to 
twenty seconds of flight time when the phugoid mode influenced 
the m.otion. After the phugoid has sufficiently subsided, the 
limit cycle motion is analagous to that previously encountered. 

To further demonstrate the character of the limit cycle 
and its independence of initial conditions, the aircraft must 
damp to a constant amplitude if the initial oscillation were 
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greater than that for the limit cycle. Figure 8 shows the 
resulting time histories for an aircraft trimmed initially 
at an angle of attack equal to 0.250 radians, well below the 
stall region. An up elevator (back stick) ramp input was 
inserted and held at a new elevator setting corresponding 
to a trim condition close to stall. Since the aircraft seeks 
the trim condition for the elevator angle prescribed, the 
final elevator angle must lie in a region analagous to the 
elevator angle for trim at stall or wherever the instability 
occurs. Because of the rapid change in elevator angle, the 
initial magnitude of the oscillation was approximately three 
times that for the limit cycle. This disturbance also in- 
troduced a large amount of phugoid motion in angle of attack 
and pitch angle as the entire limit cycle was oscillating 
with a long period movement superimposed upon the short 
period behavior. 

A ramp elevator input produced the same limit cycle. 

Figure 9 shows the time histories for the aircraft initially 
trimmed at an angle of attack equal to 0.400 radians. After 
three seconds, the elevator was deflected at a constant rate 
of one degree per second for 0.30 seconds, and then held 
constant. Again, the final elevator angle is in the region 
of the elevator angle for trim near stall. 

Figure 10 shows the aircraft trimmed at an angle of attack 
well into the stall region ( ©(,=0.4150 radians). Down ele- 
vator (forward stick) was inserted, by a step input, to an 
elevator angle near that for the stall region. The aircraft 
pitched down and the oscillation grew to the limit cycle. 
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In summary, the mathematical model used in this analysis 
has created the desired limit cycle situation. The limit 
cycle is shown to be independent of the initial conditions by 
oscillations both decaying and growing to a stable trajectory 
for the state variables. It is self-sustained as shown when 
a small impulse is inserted for 0.05 seconds and then re- 
moved. The phase plane shows a closed path trajectory, 
characteristic of a limit cycle. It must be noted that the 
actual numbers calculated are not of prime importance, rather 
the ability to simulate a known flight characteristic by 
proper aerodynamic nonlinearities was considered valuable. 

The computer program could be used for any aircraft with the 
proper data available, whether from the wind tunnel, cal- 
culated data, or actual flight test. 

IV. CONCLUSIONS 

For a multi-degree of freedom system that is nonlinear 
in its forcing function terms due to position, and has linear 
damping, it is possible to obtain a limit cycle type of 
oscillation in the system's state variables. The governing 
factor in the dynamic stability of the aircraft model was 
the slope of the lift curve. When this slope became -5.0 
or less, the aircraft was dynamically unstable in the region 
of stall. A sharp drop in lift coefficient just after 
reaching C produced this effect. This type of stall is 

characteristic of leading edge flow separation or stall. 

To produce a limit cycle motion or "rocking-chair" type 
of stall, the unstable portion of the lift curve must be 
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bounded by a stable portion on each side. With this type of 
lift curve it is possible to create a limit cycle oscillation 
in the state variables: angle of attack, pitch angle and, to 
some degree, velocity. The frequency and the characteristics 
of this oscillation were analagous to that of short period 
motion, which was predicted by analyzing the roots of the 
equations of motion. 

Since the aircraft is programmed to seek a trim condition 
at the elevator angle specified, the final elevator angle, 
after manipulation, must be close to that for trimmed flight 
in the region of maximum lift. Various types of elevator 
movement such as impulse, ramp, step, and forward stick 
produced the desired limit cycle. The resulting oscillation 
was independent of the initial conditions imposed, and dem- 
onstrated both a growing and damping oscillation to the 
same amplitude. 

The limit cycle amplitude depended directly on the length 
of the unstable portion of the lift curve. If were 

doubled, the amplitude would be doubled. This numerical 
experiment provided enough clues to encourage subsequent 
studies for analytic solutions in this area and possibly 
related topics such as post stall-gyrations, lateral and 
directional movement, and spins. 
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was held constant at -0.5730 per radian 
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FIGURE! AIRCRAFT AXIS SYSTEM 
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FIGURE 2. 

ROOT LOCUS PLOT, VARYING C, ^ 



24 




FIGURE 3a. q VERSUS CX 
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FIGURE 3c, ANGLE OF ATTACK VS. TIME 
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FIGURE 3d. LIFT COEFFICIENT VS. TIME 
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FIGURE 4. C|_ VERSUS CX 
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FIGURE 5d. PITCH ANGLE VERSUS TIME 




FIGURE 5g. velocity VERSUS TIME 
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FIGURE 6. PITCH ANGLE VS. TIME 

FOR VARYING AC|_ 
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FIGURE 9a. CK VERSUS TIME 
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FIGURE 10 b. 0 VERSUS TIME 
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THE FOLLOWING ARRAYS ARE READ INTO THE PROGRAM TO BE USED BY SUBROUTINE 
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AN ELEVATOR IMPULSE OF 1/2 DEGREE IS INSERTED FOR 1/20 SECOND 

IF( ( A.GE.2.980) -.AND. (A.LE.3.02) ) GO TO 191 
IF ( (A.GE.3.020} .AND.(A.LE.3.06) ) GO TO 195 
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RESULTANT AIRCRAFT VELOCITY ( KTS. ) VS. TIME 
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AIRCRAFT AOA VS. TIME 
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AIRCRAFT CL VS. TIME 
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AIRCRAFT CL-* (BASED ON NORMAL FLIGHT PATH ACCEL.) VS. TIME 
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THETA VS. TIME 
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